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^ I Abstract. We find a class of exact solutions to the Lighthill Whitham 

J>^' Richards Payne (LWRP) traffic flow equations. Using two consecutive lagrangian 

■ O , transformations, a linearization is achieved. Next, depending on the initial density, 

^ I we either apply (again two) Lambert functions and obtain exact formulas for the 

' dependence of the car density and velocity on x,t, or else, failing that, the same 

C/3 , result in a parametric representation. The calculation always involves two possible 

^ ^ ' factorizations of a consistency condition. Both must be considered. In physical terms, 

the lineup usually separates into two offshoots at different velocities. Each velocity soon 
(-H I becomes uniform. This outcome in many ways resembles the two soliton solution to the 

Korteweg-de Vries equation. We check general conservation requirements. Although 
traffic flow research has developed tremendously since LWRP, this calculation, being 
' exact, may open the door to solving similar problems, such as gas dynamics or water 

^ ■ flow in rivers. With this possibility in mind, we outline the procedure in some detail 

52 ! at the end. 
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^ ! 1. General history. Formulation of the model 

We have found over the years that several nonlinear, partial differential equations of 
physics, not integrable by standard methods, such as Inverse Scattering or else an 
inversion of variables, yield their secrets to lagrangian coordinate methods 
Here we will treat one such equation pair and see a combination of two 'lagrangian' 
transformations (the second one will be called quasi-lagrangian for reasons to be 
explained) and a twofold introduction of the Lambert function enable us to solve the one 
lane traffic flow problem explicitely. A further class of solutions is found in parametric 
form. 

Although we find a class of solutions that falls short of being general, common sense 
situations are well described, as well as some unexpected ones. The interesting thing, 
however, is that the much researched nonlinear equations involved, known for a long 
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time now, can in some instances be solved exactly, without recourse to approximations 
and with little or no numerics. 

In 1955, James Lighthill and his former research student, Whitham, formulated an 
equation describing single lane traffic flow, assumed congested enough to justify a fluid 
model. The theory was formulated in the second part of their classic paper on kinematic 
waves |5]. Richards independently came to the same conclusion and published in the 
following year [6]. Next Payne [7] and Whitham [8] added a second equation and replaced 
the LWR equation with standard continuity. We will call this pair LWRP. Recently the 
literature on both models has grown considerably, see for example the books by Kern 
[9] and further references p^-[T8]. 

Extensions to more than one lane, lane changing, discrete models, higher order 
effects, as well as numerical work, prevail. One of the original authors has found a Toda 
lattice like solution to the discrete version of Newell [I9] , see [20] . In a future paper, we 
will see if the methods introduced here can be applied to some of these recent extensions 
of LWRP and LWR. 



1.1. The model 

Assume a long segment of a one lane road, deprived of entries and exits, sufficiently 
congested by traffic and free of breakdowns to admit a continuous treatment, so as to 
permit us to postulate the usual equation of continuity: 

dp ^ dp du 
dt dx ^ dx 

Here p is the density of cars, the maximum of which corresponds to a rather deplorable, 
but all too familiar almost bumper to bumper situation, and u is the local velocity. 
The right hand side of the second, newtonian equation, formulated by Payne [7j and 
Whitham [8], is less obvious: 

du ^ J)u _ V{p) -u Up dp 
dt dx tq p dx' 

The first term on the right involves the mean drivers' reaction time Tq, and the next term 
models a diffusion effect depending on the drivers' awareness of conditions beyond the 
preceding car. The constant is a sort of diffusion coefficient, measuring the effect of 
the density gradient. Some recent improved models also bring in the second derivative, 
but we will not at this stage. 

In particular, u = V{p) and p both constant give a possible solution to the LWRP 
equations. Indeed, in his book Whitham expands around this equilibrium and treats the 
linear waves and possible instability so obtained |8]. However, exact nonlinear solutions 
are what we are interested in at the present. 

We specify 

V{p) = Vo - hoToy/iA^p, Vo = const. (3) 
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This often postulated form of V{p) is the only one that leads to an integrable equation, 
as far as we can see. Fortunately, the plot of Q{p) = pV{p), which people measure, is 
always convex and parabola-like, see Whitham's figure 3.1. Thus, so far luck is with us. 
We introduce dimensionless variables by replacing 

t tTo, {u,Vo) ^ {u,Vo)y/iA^, X x^/z^To, p(/ioro)"^ (4) 

This leaves the continuity equation unchanged, and the newtonian equation takes the 
form: 

du du ^ 1 dp . 

dt ^ dx ^ ^ ^ p dx 

2. Introducing lagrangian coordinates 

The non-linearity on the left hand side of equations ([1]) and ([5]) can be eliminated by 
introducing lagrangian coordinates: ^(x, t), the initial position (at t = 0) of a fluid 
element which at time t was at x, and time t. In this description, the independent 
variable x becomes a function of ^ and t, as are the fluid parameters p(^, t) = p{x{^, t),t) 
and u{^,t) = u{x{^,t),t). 

Here and in what follows we adopt the convention that a superposition of two 
functions which introduces a new variable is denoted by the same symbol as the original 
function, but of the new variable. Denoting by / either p or u, the basic transformation 
between eulerian coordinates x, t and lagrangian ones ^, t can be written as 

,^ r f,,'^A,' f,,^ dfix,t) , dfix,t) 

x{^,t) = ^+ u{^,t)dt, —=u{^,t), — - — = — - — + u — . (6) 



dt ' ' dt dt dx 

We denote by s{^) the number of cars between the last one at ,^ = C^^ain and that at ^: 

^(0 = r PoiC) dC, Po(0 =p{x = U = 0) (7) 

^ ^min 

where Po(0 is the initial mass density distribution. 

Here and in what follows, the subscript always refers to t = 0. We will also use 
the superscript to refer to ^ = 0. 

There are no gaps in the line of traffic considered and therefore Po(0 positive. 
Hence s(^) is an increasing function starting at s(^inin) = 0, and one can introduce a 
uniquely defined inverse function ^(s). The initial position of a fiuid element can be 
specified by either ^ or s. 

If a small initial interval d^ at t = becomes dx at time t, mass conservation 
requires: 

ds = po(Od^ = p{x,t)dx. (8) 
This leads to a mass conservation equation in lagrangian variables: 

dx{s,t) ^ 1 . , 

ds p{s,t)' ^ ' 
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and to a useful operator identity 

= = (10) 

p{x,t)dx PoiO 9^ 9^ 

Integrating (jH]) over s' from s{C, = 0) to s, we obtain the continuity equation in 
integral form: 

ds' 



X{s,t) = x{s,t) 



11 



sO Pis',t) 

This indicates that if we know the car density in lagrangian coordinates p{s,t), we can 
determine the evolving shape of the line of traffic, where the distance X is measured 
from the ^ = car. 

The analog of the continuity equation ([T]) is obtained by differentiating by t. 
Using the middle part of we obtain 

d^jj{s,t) du . 1 

ds ' 



(12) 



dt OS p 
The newtonian equation in lagrangian coordinates is obtained from and flTU]) : 



du{s,t) dp 



(13) 



Equation ( IT3l) is linear and can be solved to express u{s,t) in terms of p. Again, using 
the middle part of ([6]) we can also calculate x{s,t): 

ft 



u{s,t) = e ' 
N{s,t) = Vo 



N{s,t')e'' dt' + u{s,0) 



P + 



dp 
ds 



x{s,t) = ^(s) + / u{s,t')dt' 
Jo 

= ^{s)+u{s,0)-u{s,t)+ I N{s,t')dt', 

Jo 

where the function u{s, 0) will be determined later. 



(14) 
(15) 

(16) 



3. Finding the fluid density 

Differentiating the newtonian equation (fT3l) by s, and using continuity (fT2|) . we obtain 
one equation for ip: 

d'^tfj d / 1 dtlj\ dtp d I Q 
df^ dsyip'^dsJ dt ds 
This equation can be factorized in two possible ways, I and II: 



17) 



I • f ^ + M f ^ + l + ip] - 

\dt ds thJ \ dt lb ds J 



and 



dt ds TpJ \ dt ijj ds 
d d 1\ {dib 1 dip 



Exact solutions to the Lighthill Whitham Richards Payne traffic flow equations 5 



We will find that the second factor in (1181) best yields solutions such that X > 0, 
whereas that in f|T9|) rules X < 0, where X is always the distance from the car that 
started at x = 0. Different pairing would lead to trouble. 

If we denote equation (fT7|) as Oip = and the decompositions by I and II, the 
following symmetries hold: 

= OlOl = 0'^0'l olii;) = i-iy^'Ofi-i;), z = l,2. (20) 

In what follows, we will find solutions for which the second factor in one of equations 
( |T8|) . ( fT9i) vanishes, leaving a more general treatment to a possible later paper. We follow 
motion from left to right. Factorization also means that we can only introduce the initial 
value of the density (or tp). The initial velocity u{s,t = 0) will then follow except for a 
universal constant. We will have more to say about this later on. 

The non-linearities in f|T8|) and f|T9l) (second factors) can be eliminated if one 
transforms the variables s,t to ri,t in a way similar to the lagrangian transformation 
IQ, though without the usual physical interpretation: 

s[r],t)-7]T ^^^^^,y ^ ^^p{rj,t)' dt dt ^^ds- ^ ^ 

Solving the resulting linear equation 

ot 

we obtain, in view of the fact that s and rj are identical at t = 0, 

^(r/, t) = Tl + ^AMv) ± 1], Mv) = ^{s = r],0). (22) 
For this iplrjjt) we have 
"* df 

— = TH^'Hv,t)/Mv)], (23) 

and finally, back to p = l/ip and using fl2T]) . 

s = r] + \n{lTPo{v)Mt)), po{v) = Po{s = v) , A{t) = - 1. (24) 

In this relation, defining s in terms of rj and t, Po{t]) is defined by (IT!) but is expressed in 
terms s, where one has to rename s to rj. Exactly the same procedure applies to ipoiv) 
given by ([22]). 

Using fl22|) . we can express p in terms of rj and t: 

Piv,t) = , , . , (25) 

=Fl + e ^[l/po{v> ± IJ 

which tends to po(s) as t — )■ 0. 

We are now in a position to determine the function u{s, 0) needed in equations dHj) 
and (fT6l) . Differentiate (fT3l) by s and then subtract both sides of (fT7|) from the result to 
obtain 

d \ /dip du\ 



Exact solutions to the Lighthill Whitham Richards Payne traffic flow equations 6 



Solved by 



Therefore, if f{s) = 0, equation (|T2l) will be valid for all time. All we require is 



'dt 



ds 



0, 



I.e. 



du{s, 0) dipo 



t=o 



ds 



dt 



(27) 



(28) 



This result, along with either ( !T8|) or (fT9|) . leads to 

du{s,0) ^ _^]_di^ T^-i^o 
ds ipo ds 

Integrating over s' from s° to s and transforming the result to C,, we end up with 

Po(0- 



s' = si^ = 0), 



Po(e = 0), (29) 



where Mq = u(^ = 0, 0) > is arbitrary. 

The last task is to determine u{s,t), x{s,t), and X{s,t), given by flTl]) . ([TB|1 . and 
( ITTj) . in terms of r^. Using ( p5|) . (12^ and ( !T5|) we find the integrand iV: 



iV(^5,t) = K)- 



9p 



Vr 



P + 



dp/drj 



±1 



ds/drj 
Po{v)+Poiv) 



which tends to Vq ± 1 as t -t- oo. Here rj = ri{s,t) must be found as a solution of 
the transcendental equation (^^, and the integrals flT^ and ( IT6|) must be calculated 
numerically. On the other hand, the integral (ITT]) can be calculated analytically: 



X(s,t) = x(s,t) - 



x(s^^) 



p(v,t) 



dr/' 



77) - ^(s = 77°) =F A(t) ?7 - 77° + In 



Po(^) 
Po(r/°) 



(31) 



where rj = ri{s,t) and 77° = 77(5°, t) are defined implicitly by ( l24l) . 

Notice that the time evolution of an assumed initial density profile po(0 can only 
be determined if the solution r]{s,t) of equation (IMl) is a continuous function of t. This 
will be true either in case I, if Po(0 decreases from its initial value a = po(^ = 0) as 
increases (i.e., for C, > 0), or in case II, if Po(0 grows as C, increases from —00 to 0. 

In either case, f]{s,t) starts its time evolution with ri{s,0) = s. Then in case I, it 
grows to rjranxis) = ri{s,t — > 00), and in case II, it falls to ri^in{s) = ri{s,t 00). 

If po(0 has a maximum at ^ = either for ,^ > or ^ < 0, then for any ^ within the 
interval (0,^m), ri{s{^),t) as a function of t would have a jump, which is unacceptable. 

Another obvious requirement on the initial density profile po(0 is its integrability 
over the interval (—00, 0) or (0, 00). 

In practice, the general theory given in the last two sections is only useful if the 
integral (171) can be determined analytically. Examples will be given in the following 
sections. 
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4. Two exponential profiles of the initial fluid density 

We will see that two exponential profiles of the initial fluid density 

Po(0 = «exp(-AO, e>0, i.e. Un = 0, (32) 
Po(0 = «exp(AO, ^ < 0, i.e. ^min = -oo, (33) 

play a special role here, as in their case it is possible to eliminate the auxiliary variable 
7], and even find the fluid density p in terms of X and t. This is because the argument 
of the logarithm in equation ( 12^ defining 77(5, t) is a linear function of 77. Therefore the 
solution ri{s,t) of (12^ can be given in terms of the Lambert function W{x) defined by 

Wexp{W) = X, equivalent to W + \iaW = Inx, -1/e < x < 00. (34) 

The W function is named after Johann Heinrich Lambert, Leonhard Euler's young 
protege and an important mathematician in his own right (1728-1777), see [21] for a 
summary. For negative x, W has two negative branches and the upper one (continuous 
at a; = 0) along with the main branch of the logarithms should be chosen in what follows 
{Inx = In |x| + Itt when x < 0). 

The solution Y of the transcendental equation 

Y + \n{aY + b) = c (35) 

is given by 

r = iy[expfc+ ^ -Ina)] - -. (36) 
L V " /J a 

Using equation ([7j) we first find 

ro ( for ([32D, 

,o = ,(e = 0)= / Po{ade={ a , _ (37) 

Cm i n 



for 



and then calculate 



^(l-exp(-AO) form. 

^exp(AO for ([33]). 

A 



(38) 



The inverse functions are given by 



1 , , I -ilnfl - for 

■In— for ([33]). 

A a 

Using this formula we can transform the initial conditions ( l32l) and ( l33l) given above in 
X, if: to s, t: 




, ^ ... for (E2), 
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We now look for solutions to equations (ITS]) and (HM that recreate the above initial 
conditions as t tends to zero. 

Replacing s by in f HOl) and using the po{f]) so obtained in ( 12^ . we can write the 
result as 



r/ + ln|AA(t)r7+ 1 + A(t)(Ta - As°) |. 



(41) 



Indeed the argument of the logarithm is a linear function of t], and we can use equations 
( 135|) and ( 136|) to obtain 



T]{s,t) = W{s,t) ± - + S° - 

A 

w{s,t) = M/^jexp p{s,t) |, 



1 



\A{t)^ 



A{t) 



(42) 



p{s,t) 



^ A Av4(t) 



+ In 



XA{t)' 



Inserting this f]{s,t) into (!25|) we obtain, in view of f HOj) 



p{s,t) 



(43) 



/l(t) lXA{t)w{s,t) 

We recover the initial condition (140|) when t — )■ 0+, and so 1/A(t) — ?■ oo, by invoking 
the large argument approximation for W: W{x) ~ ln(x) — lnln(a;), x ^ 1, 



A 



- 1 



a =F A(s — s*- 



1 + A{X{s - s^) T a)) 



lXA{s - 3° T f + ^) 
= aTX{s~s°) + 0{A). (44) 

Using (143|) and ( |29|) we can determine A^(^, t) and m(^, 0) needed in equations f lT4|) - 
( !T6ll . where s = s(^) is given by ( l38ll : 



iV(e,t) = K)±l± 



1 



1 



1 



M(e,0) =no + (A-l)e 



A(l -e-*)(l + w(s,t)) 

as t — 7- oo, 
a 



X 



exp{TXO 



(45) 
(46) 



where in equation (H5|) we used the fact that for large t, w{s,t) ^ — ln[A(e* — 1)] which 
tends to minus infinity as t — )■ oo. Calculating the integrals in ( |T4|) and ( |T6l) numerically, 
we can find characteristics x(^, t) parametrized by the initial fluid element position ^, 
shown in figures 2 and 4. 

We will now try to express p directly in terms of X by using ( ITTl) . Inserting p(s, t) 
given by f l43|) into f lTT]) and changing the integration variable from s to we obtain 

w{s',t) 



XX = ±X^A\t)e-' [ ds' 

JsO 



XA{t)w{s',t) 



±XA{t)e-' 

±e-* 



w{s,t) 



dw' 



1 + 



w{sO,t) 

XA{t)[w{s,t) - w;(s°,t)] 

XA{t)w{s,t) - 1 



1 + XAjt) 
XA{t)w' - 1 



[1 + XA{t)] In 



AA(t)u;(sO,t) - 1 



(47) 
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which can be written as 

Y + ln([l + XAm - l) - ^'^^VAl(!r^'°'^ + ln(AA(t).(.o, t) - 1 



+ \A{t) 
where 

\A{t)w{s,t) 
l + \A{t) ■ 

Again using (!35l) and (l36ll to solve for F, we find Ay4(t)w(s, t) as a function of X and t. 
Inserting it into fH3|) we end up with 



[l + \A{t)]W{et>) 

p{X,t) = {^{l-e-^){l + [l + \A{t)]W{et>)y " ' (48) 

0, X<0, 
where 

TAe^X + AA(t)M.o,t)-l , , AA(t)M.°,t)-l 

= rTX4(t^ + 1 + AA(t) • ^'9) 

Note that equation (1471) could also be obtained from ( 13T1) . where is given by (l42l) . 

in view of ( 139|) . and po is given by (HOjl with s = rj. 

As t — )■ we recover the initial condition (1321) or ( l33l) by once again invoking the 
large argument approximation W{x) ~ ln(a;) — lnln(a;). The calculation is similar to 
that in s, t above. 

Jf^.l. Exponentially decreasing initial fluid density 

As initial condition at t = we first take a pulse described by ( |32|) and subsequently 
expected to move to the right. To the far right cars are thin on the ground. We expect 
them to move freely and make it possible for the congestion near the discontinuity to 
spread out. This common sense expectation will test our result. 

The density profile p(X, t), defined by equation (148|) with upper sign is shown in 
figure 1. 

As t — 7- oo this profile tends to 

r [l-exp(-a/A)]exp(-X) 

J l-[l-exp(-a/A)]exp(-X)' " ' (50) 

[0, X < 0. 

That is because AW ~ —[1 — exp(— a/A)] exp(— X), as small argument W functions can 
be approximated by their argument, W{x) ~ x, regardless of sign. The initial condition 
only appears as the ratio a/ A. A whole class of initial profiles ends up developing 
identically, as long as this ratio is the same (figure 1). 

In figure 2 we present characteristics a;(^, t) obtained by numerically calculating the 
integrals in (fl^ and (|T6|) . No crossing of characteristic curves is observed, a pleasing 
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Figure 1. Two density profiles for various times as found from our solution (I). Here 
a = 1, A 2 in the first case, and a = 4, A = 8 in the second one. Nevertheless, the 
emerging profiles are seen to be identical after a while. The value of a for each surface 
can be seen as equal to p(0, 0). 




Figure 2. Characteristics x{^,t) as functions of t for ^ ranging from to 0.2 in steps 
of 0.05, a = 4, A = 8, ito = 1, and Vq = 2. 

vindication of the LWPR model. In particular, the ^ = line in figure 2 gives us the 
trajectory of the discontinuity point, defining X. 

Equation ( H5ll with upper sign implies a uniform motion at the final, large t stage, 
at which u = Vq + 1. This can also be seen from equation ([2]). This emerging solution is 
identical to that given by Whitham |8j as a special case of a 'continuous shock structure', 
Whitham's equation (3.16), (when his A = 0,U — v = 0). Of course, continuity is lost at 
s = 0. Now we know exactly how to set up initial conditions so as to obtain this profile. 
Alternatively, we could set up this profile from the beginning. There are obviously 
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2.0 




2 PCS.t) 



3 



4 



-2.0 



X 



Figure 3. Two density profiles for various times as found from our solution (II). Here 
a = 1, A = 2 in the first case, and a = 4, A = 8 in the second one. Nevertheless, the 
emerging profiles are seen to be identical after a while. The value of a for each surface 
can be seen as equal to p{0,Q). Notice that the final amplitude of p is smaller here 
than in figure 1. 



fewer possible final states than initial conditions, a state of affairs often encountered in 
non-linear problems, see e.g. [221 [23], [2] for further references. 

We can now check to see if the integral of p is conserved in the t — ?■ oo limit. In 
view of dSO]) 



^0 

The integral over all X of p{X,t — )■ oo) is indeed a/A. Our exact solution converges to 
a stationary mass conserving travelling wave. 

4.2. Exponentially increasing initial fluid density 

We now find a similar solution, but to case II. Our initial condition fl33|) is now limited 
to a; < 0. 

This is the mirror image of the previous lineup, with the most congested traffic 
facing an empty road. However, motion will still be from left to right. 

The density profile p{X,t), defined by equation f l48p with lower sign is shown in 
figure 3. As t — )■ 00 this profile tends to 




which implies 





X < 0, 



(51) 
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Figure 4. Characteristics x{£^,t) as functions of t for ^ ranging from to —0.025 in 
steps of -0.005, a = 1, A = 2, uo = 0.05, and Vo = 10. 



The end result is very similar to the previous one, but not identical. Mass 
conservation follows from the same calculation. However, we are solving a different 
equation and the asymptotic, uniform velocity is Vq — 1 as opposed to the Vb + 1 of case 
I, see (H5l) with lower sign. This can also be seen from figure 4, where characteristics 
x{^, t) are presented. 

If we match the initial velocities at the centre, x = 0, then after a while, two 
cavalcades emerge, a faster one at velocity u = Vq + 1, and a slower one at Vq — 1 
(cheaper cars?) see figure 4. Our situation reminds us of solutions to the wave equation. 
Here {cdx — dt)ip = would correspond to I, {cdx + dt)ip = to II. Coexistence of the 
two was not surprising in the linear world, now it is somewhat unexpected. Although 
we are combining two solutions following from different factorizations of the governing 
equation, we should remember that we are nevertheless dealing with one exact solution, 
unique to the initial profile and velocities. 



5. The initial density profiles that can be treated parametrically 

In this section we present a few initial density profiles satisfying the applicability 
conditions of our theory as formulated at the end of section 3, see figure 5. 
Detailed calculations will be performed for a pair of cases: 

Po(0 = ^^^^^^ = « [1 - tanh2(A0], (52) 

where either < < oo in case I, or — oo < < in case II. 

The fact that the derivative dpo(0/d^ vanishes at ^ = 0, in contrast to the 

exponential profiles (132|) and ( !33|) . will have a consequence on the time evolution in 
case I, see figure 6. 



Exact solutions to the Lighthill Whitham Richards Payne traffic flow equations 13 



Normalized density 
l.O-t 




Normalized position 



Figure 5. Normalized density profiles po — po/a versus normalized position ^ = 
for poiO given by (a): (|32]) and (b): ^ for 6 = 1, r = 3, (c): ([52]), and (d): 

(El. 



The remaining profiles will have a power behaviour at infinity, po(0 ~^ (=^0 ^ &s 
— )■ oo, where r is a real number greater than unity for integrability: 

a 



PoiO 



and 



Po(0 = « 



1 + (A0'' 



6 > 0, 



r > 1, 



(53) 



(54) 



(±Ae + &)'■' 

where the upper sign refers to case I, ^ > 0, and the lower one to case II, ^ < 0. 

By analogy to the exponential profiles fl32p and fl33p . each pair of symmetric cases 
can be treated in a single calculation. For po(0 given by fl52|) we first find 



5° = s(e = 0) 



(0 for ^ > 0, 

(min for ^ < 0, 



and then calculate 

s{0 = ^° + / Po(^') = ^' + J tanh(AO. 
The inverse functions are given by 

{1 + 
2A 1 - As/a 
1 Xs/a 
2A 2 -Xs/a 

Using tanh(A^) calculated from (!56l) in (!52l) we obtain 

2i I a[l - (As/a)2] 



(55) 



(56) 



Po[s) = a 



1 



A(s-s°)/a 



As(2 - As/a) 



for ^ > 0, 
for ^ < 0. 

for ^ > 0, 
for ^ < 0. 



(57) 



(5^ 
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Replacing here s by 77 and using the po(^) so obtained in and along with (15^ we 
find equations defining ri{^,t) and the integrand N{r],t) needed in equations f|T4|) - f|T6|) : 

for ^ > 0, (59) 
for ^ < 0, (60) 



J tanh(AO = r/ + ln(^l - a[l - {Xr] / a)^]A{t) 

tanh(AO] = + ln(^l + A r/(2 - Xr]/a)A{t) 
Tl-fiv) 



N{r],t) = Vo±l 



iTf{v)A{ty 



where 



fiv) 



-r]{7] + 2)Xya 



X 



-T]^X/a + 2r]{l - A/a) + 2 



for ^ > 0, 
for ^ < 0. 



(61) 



(62) 



In a similar way we can determine X{ri,t) by using flHT]) along with flFTl) and 
with s = r]-. 

"* " ' l^- l-Xr]/a f.^,. 1\, l + A?7/a 



X = — 



X 



XA{t) 



In 



+ XA{t)iv-v') 



1 - Xriya 
for ^ > 0, 



XA{t) 



In 



1 + Xrj^a 



(63) 



+ AA(t)(r/-r^°) 



for ^ < 0. 



(64) 



Using here ?7(^, t) defined implicitly by either of equations fl60|) and in p{ri, t) given by 
( !25|) . we obtain p(X, t) in parametric form: p(^, t) and X(^, t). This form is appropriate 
for making use of the ParametricPlotSD command of Mathematica. The results are 
shown in figures 6 and 7. They resemble those shown in figures 1 and 2 except for the 
time evolution of the discontinuity at X = 0, p{X = 0,t) for case I (^ > 0), shown in 
figure 6. 

The characteristics x{^, t) can be found from equations (n31) -( fT6ll by numerical 
integration, where the integrand N{^,t) is defined by (!6T|) and either of equations (160|1 . 
and 

a 



u{^, 0) = so-^t[j tanh(AO - 2 In cosh(AO } 



(65) 

see equations (!29| (!52|) and (156|1 . The results, depending on two parameters Vq and uq, 
are shown in figure 8. 

A characteristic feature of the plots representing the density given in parametric 
form, p{^,t) and X{^,t), is that the mesh lines correspond to ,^ = const, and t = const, 
see figures 6 and 7. For the density given explicitly, p{X,t), they correspond to 
X = const, and t = const, see figures 1 and 2. The mesh lines ^ = const are particularly 
useful. Each point on a. ^ = const mesh line gives us the information of both the actual 
position X and the associated density at time t, for the car that started from X = ^ 
at t = 0. This information is given in the frame moving along with the discontinuity at 
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2.0 




2.0 



Figure 6. The fluid density represented parametrically as found from our solution (I). 
Here a = 1, A = 2. Tire mesh lines correspond to t = const or ^ = 0, 0.25, 0.5, 0.75, . . .. 
Each value of ^ is equal to X at i = 0. 




Figure 7. The fluid density as in flgure 6 but as found from our solution (II). Here 
again a = 1, A = 2, and ^ = 0, —0.25, —0.5, —0.75, . . ., see X at t — 0. 

^ = 0. The motion of these frames in turn is described by the characteristics labelled 
with ^ = in figure 8. 

Adding cases I and II, we have a solution such that the initial congestion splits in 
the middle, resulting once again in a slower cavalcade following a faster one, see figure 
8. This is rather like a two soliton solution of the Korteweg-de Vries equation, see e.g. 

m- 

6. Summary 

Our solutions augment those found for simpler, single equation nonlinear models, e.g. 
Burgers, see [8J. Our exact solutions converge to single or double stationary travelling 
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6 - 

^=0.5 




Figure 8. Characteristics x{£,,t) as functions of t for ^ = 0, 0.25, 05, and ^ = 
0, -0.25, -05, a = 1, A = 2, Mo = 7, and Vq = 10. 

wave structures after a few tq (figure 1). This steady state convergence is more or less 
wliat one would expect, at least in the single case. The LWRP equations have several 
families of solutions, to which we hereby add. Generally we must use common sense to 
home in on the relevant ones. However, exact solutions help and absolve us from the 
need to guess. As in police work, connecting initial conditions with a present situation 
requires painstakingly following through the history. Exact solutions exempt us of this 
task. 

We introduce a lagrangian transformation in x,t. The variables are now C,,t, ^ 
being the position at t = 0. It labels the car rather than its actual position. The 
continuity equation is now linear, simplifying further calculations. It may so happen 
that the newtonian equation can yield a condition in just one dependent variable. If this 
condition can be factorized, there is a possibility of finding exact solutions by finding 
solutions to this factorized equation, being of lower order. In our case, two distinct 
factorizations were possible, both were solved, and both solutions were needed to cover 
all X. In fact, each factorization was assigned to a specific half axis, the same for all 
cases considered. This situation, demanding merging two solutions to cover the whole 
domain, is unusual in a nonlinear problem. Another requirement was that the original 
profile have only one maximum at x = 0. 

We had to introduce a second, quasi lagrangian transformation to linearise. In 
general, we obtain the solution in parametric form. However, if we are lucky, we can 
solve completely, that is get rid of the parameter, by using Lambert functions. We give 
an example of this. 

It should be stressed that a complete solution is only possible if we combine our 
two factorized equations, I and II. This is similar to the situation for solutions to the 
wave equation 
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such that solutions of {cd^ — dt)ip = and {cdx + dt)ip — coexist. However, this 
coexistence is somewhat surprising for a nonhnear system. 

An interesting question is: how wide a class of non-linear problems can be so solved? 
Perhaps one possibility is furnished by the somewhat similar gas dynamics and shallow 
water equations? In our case, exponential initial conditions alone seem to promise that 
the Lambert function will make an appearance. Perhaps they should be used for similar 
problems. 

Another interesting question is how general is the requirement of merging solutions 

I and II when two factorizations of the consistency condition exist. 
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